/*******************************************************************************

This code file produces Figure A9, "Event Studies of GEA and NPP Reforms."

*******************************************************************************/

*** Manage settings

	run "$dir/code/modules/settings.do"
	
*** Run LPMs, store coefficients in matrix

	use "$data/clean/cleaned_data.dta", clear
	
	matrix results = J(13,9,.)
	
	** GEA
	
		* Tax incentive 
		reg dtaxrate_onsite 1.gea_expansion#i.yearpermit i.yearpermit $controls_lot $controls_block i.borough i.yearpermit, cl(nta)
		
		forvalues i = 1/13 {
			
			local y = `i'+2002
			
			lincom 1.gea_expansion#`y'.yearpermit - 1.gea_expansion#2007.yearpermit
			
			local tmp1 = r(estimate)
			local tmp2 = r(se)
			
			mat results[`i',1] = `y'
			mat results[`i',2] = `tmp1'
			mat results[`i',3] = `tmp2'
			
		}
		
		* Has inclusionary units
		reg inclusionary_onsite 1.gea_expansion#i.yearpermit i.yearpermit $controls_lot $controls_block i.borough i.yearpermit, cl(nta)
		
		forvalues i = 1/13 {
			
			local y = `i'+2002
			
			lincom 1.gea_expansion#`y'.yearpermit - 1.gea_expansion#2007.yearpermit
			
			local tmp1 = r(estimate)
			local tmp2 = r(se)
			
			mat results[`i',4] = `tmp1'
			mat results[`i',5] = `tmp2'
			
		}
	
	** NPP
	
		* Comparison group is Census tracts that overlap with an NPP area
		bys ct2010: egen sh_npp = mean(npp)
		keep if sh_npp > 0
	
		* Tax incentive 
		reg dtaxrate_onsite 1.npp#i.yearpermit i.yearpermit  $controls_lot $controls_block i.borough i.yearpermit, cl(nta)
		
		forvalues i = 1/13 {
			
			local y = `i'+2002
			
			lincom 1.npp#`y'.yearpermit - 1.npp#2007.yearpermit
			
			local tmp1 = r(estimate)
			local tmp2 = r(se)
			
			mat results[`i',6] = `tmp1'
			mat results[`i',7] = `tmp2'
			
		}
		
		* Has inclusionary units
		reg inclusionary_onsite 1.npp#i.yearpermit i.yearpermit $controls_lot $controls_block i.borough i.yearpermit, cl(nta)
		
		forvalues i = 1/13 {
			
			local y = `i'+2002
			
			lincom 1.npp#`y'.yearpermit - 1.npp#2007.yearpermit
			
			local tmp1 = r(estimate)
			local tmp2 = r(se)
			
			mat results[`i',8] = `tmp1'
			mat results[`i',9] = `tmp2'
			
		}
	
*** Move results to data

	svmat results
	keep results*	
	rename (results1 results2 results3 results4 results5 results6 results7 results8 results9) (year dtaxrate_coef_gea dtaxrate_stderr_gea inclusionary_coef_gea inclusionary_stderr_gea dtaxrate_coef_npp dtaxrate_stderr_npp inclusionary_coef_npp inclusionary_stderr_npp)
	
	foreach s in gea npp {
		
		replace dtaxrate_stderr_`s' = . if year == 2007
		replace inclusionary_stderr_`s' = . if year == 2007
			
		gen dtaxrate_lo_`s' = dtaxrate_coef_`s' - 1.96*dtaxrate_stderr_`s'
		gen dtaxrate_hi_`s' = dtaxrate_coef_`s' + 1.96*dtaxrate_stderr_`s'
		
		gen inclusionary_lo_`s' = inclusionary_coef_`s' - 1.96*inclusionary_stderr_`s'
		gen inclusionary_hi_`s' = inclusionary_coef_`s' + 1.96*inclusionary_stderr_`s'
		
	}
		
*** Make graphs 

	** GEA

	tw rcap dtaxrate_lo_gea dtaxrate_hi_gea year, lcolor(gs9) || line dtaxrate_coef_gea year, lcolor(navy) ///
		graphregion(color(white)) xtitle("Year of Permit") legend(order(2 "Estimate" 1 "95% CI")) ///
		ytitle("421-a Tax Incentive") name(gea_dtaxrate, replace) xscale(range(2003 2015)) xlabel(2003(2)2015) legend(size(small)) yscale(range(-0.1 0.5)) ylabel(-0.1(0.1)0.5)
		
	tw rcap inclusionary_lo_gea inclusionary_hi_gea year, lcolor(gs9) || line inclusionary_coef_gea year, lcolor(navy) ///
		graphregion(color(white)) xtitle("Year of Permit") legend(order(2 "Estimate" 1 "95% CI")) ///
		ytitle("Has Onsite Inclusionary Units") name(gea_inclusionary, replace) xscale(range(2003 2015)) xlabel(2003(2)2015) legend(size(small)) ///
		yscale(range(-0.05 0.7)) ylabel(0(0.2)0.6)

	** NPP

	tw rcap dtaxrate_lo_npp dtaxrate_hi_npp year, lcolor(gs9) || line dtaxrate_coef_npp year, lcolor(navy) ///
		graphregion(color(white)) xtitle("Year of Permit") legend(order(2 "Estimate" 1 "95% CI")) ///
		ytitle("421-a Tax Incentive") name(npp_dtaxrate, replace) xscale(range(2003 2015)) xlabel(2003(2)2015) legend(size(small))
		
	tw rcap inclusionary_lo_npp inclusionary_hi_npp year, lcolor(gs9) || line inclusionary_coef_npp year, lcolor(navy) ///
		graphregion(color(white)) xtitle("Year of Permit") legend(order(2 "Estimate" 1 "95% CI")) ///
		ytitle("Has Onsite Inclusionary Units") name(npp_inclusionary, replace) xscale(range(2003 2015)) xlabel(2003(2)2015) legend(size(small))
	
	** Combine
	
	grc1leg gea_dtaxrate gea_inclusionary npp_dtaxrate npp_inclusionary, cols(2) rows(2) ///
		graphregion(color(white))
		
	graph export "$figs/event_study_gea_npp.pdf", replace
	graph export "$figs_overleaf/event_study_gea_npp.pdf", replace
